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I Milagro is a water Cherenkov extensive air shower array that continuously 

^ I monitors the entire overhead sky in the TeV energy band. The results from 
an analysis of ~3 years of data (December 2000 through November 2003) are 

Tjj- ' presented. The data has been searched for steady point sources of TeV gamma 
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rays between declinations of 1.1 degrees and 80 degrees. Two sources are detected, 
the Crab Nebula and the active galaxy Mrk 421. For the remainder of the 
Northern hemisphere we set 95% C.L. upper limits between 275 and 600 mCrab 
(4.8-10.5x10^^^ cm~^ s^^) above 1 TeV for source declinations between 5 degrees 
and 70 degrees. Since the sensitivity of Milagro depends upon the spectrum of 
the source at the top of the atmosphere, the dependence of the limits on the 
spectrum of a candidate source is presented. Because high-energy gamma rays 
from extragalactic sources are absorbed by interactions with the extragalactic 
background light the dependence of the flux hmits on the redshift of a candidate 
source are given. The upper limits presented here are over an order of magnitude 
more stringent than previously published limits from TeV gamma-ray all-sky 
surveys. 

Subject headings: gamma rays: observations — surveys — galaxies: active 

1. Introduction 

Sources of very-high- energy (VHE, >100 GeV) gamma rays arc observed to be non- 
thermal in nature and arc typically the sites of particle acceleration. This acceleration is 
thought to occur in astrophysical shocks such as those believed to exist in plerions (De 
Jager and Harding 1992), supernova remnants (Berezhko and V51k 2000), active galactic 
nuclei (Blandford and Ostriker 1978), and galaxy clusters (Loeb and Waxman 2000) (among 
other sources). These shocks may accelerate protons or electrons, both of which lead to the 
emission of gamma rays. Since gamma rays are unaffected by the magnetic fields that pervade 
the Galaxy and the Universe, they can be used to pinpoint the sites of particle acceleration. 
In addition to these "classical" astrophysical sources of VHE gamma rays, other more exotic 
objects such as primordial black holes, topological defects, and the decay of relic particles 
from the big bang may also emit VHE gamma rays. Perhaps of most interest is the possible 
existence of a new type of source that has yet to be postulated. A comprehensive survey 
of the sky sensitive to emission at all time scales is necessary to detect the many possible 
sources. The analysis presented in this paper is part of an ongoing effort by the Milagro 
collaboration to search the entire northern hemisphere for such objects. The search for short 
bursts of TeV gamma rays has been addressed in a previous paper (Atkins et al. 2004). The 
analysis presented here deals specifically with steady point sources of TeV gamma rays. 

It has been over a decade since the discovery of the first source of VHE gamma rays, the 
Crab Nebula (Weekes et al. 1989). Since that time there have been seven other confirmed 
sources of TeV gamma rays (Horan and Weekes 2003), six of which he in the northern 
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hemisphere. With the exception of the Crab Nebula and the supernova remnant PKS 1706- 
44, these objects are all active galactic nuclei of the blazar class (Horan and Weekes 2003). 
All of these objects have been discovered by atmospheric Cherenkov telescopes searching for 
counterparts to sources discovered at lower energies. In contrast, the EGRET instrument 
detected over 270 objects emitting high-energy gamma rays above 100 McV (Hartman et al. 
1999). One hundred and seventy of these objects are not identified at other wavelengths. 
The VHE regime is a natural energy band to search for counterparts to these objects. The 
small field-of-view and low duty factor of the atmospheric Cherenkov telescopes (ACTs) 
make comprehensive sky surveys difficult to perform. As a result, very few comprehensive 
surveys of the VHE sky have been performed to date. The first VHE survey was performed 
by a non-imaging ACT (Helmken et al. 1979; Weekes et al. 1979) in 1979. The Milagrito 
instrument (a prototype to Milagro, with no background rejection capability and a higher 
energy threshold) also performed a survey of the northern hemisphere (Atkins et al. 2001) 
and set limits of <3 Crab from any point source in the northern hemisphere. More recently 
flux limits of 4-9 times that of the Crab Nebula (above 15 TeV) have been published by the 
AIROBICC collaboration (Aharonian et al. 2002). (1 Crab is equivalent to an integral flux 
above 1 TeV of F(>1 TeV)= 1.75x10"^^ cm~^ s~^.) The hmits presented here are over an 
order of magnitude more stringent than these previous surveys. 

Since many of the conflrmed sources of VHE gamma rays are extragalactic the hm- 
its must account for the absorption of TeV gamma rays by interactions with the extra- 
galactic background hght (EBL) (Primack et al. 2000; Stecker and De Jager 2002; Kneiske, 
Mannheim, and Hartmann 2002). The EBL is comprised of visible radiation emitted by 
stars and infrared radiation emitted by dust due to reprocessed starlight. Since direct mea- 
surements of the intensity and spectrum of the EBL arc problematic due to the foreground 
light from our galaxy, a model is used to determine the effect of the EBL on the observed 
spectrum at earth from a distant source. 

Before employing the results presented here the limitations of this survey need to be 
understood. First, the limits presented only apply to point sources, not to extended objects, 
such as the galactic plane. Second, they only apply to the average VHE emission during the 
time period over which the data was obtained. 

2. The Milagro Observatory 

The Milagro gamma-ray observatory has 723 photomultiplicr tubes (PMTs) submerged 
in a 24-million-liter water reservoir. The detector is located at the Fenton Hill site of Los 
Alamos National Laboratory, about 35 miles west of Los Alamos, NM, at an altitude of 2630 
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m asl (750 g/cm^). The reservoir measures 80m x 60m x 8m (depth) and is covered by 
a hght-tight barrier. The PMTs are arranged in two layers, each on a 2.8m x 2.8m grid. 
The top layer of 450 PMTs (under 1.4 meters of water) is used primarily to reconstruct the 
direction of the air shower. By measuring the relative arrival time of the air shower across 
the array the direction of the primary cosmic ray can be reconstructed with an accuracy of 
roughly 0.75°. The bottom layer of 273 PMTs (under 6 meters of water) is used primarily 
to discriminate between gamma-ray-initiated air showers and hadronic air showers. 

The discrimination of the cosmic ray background is described in detail in (Atkins et al. 
2003). The background rejection uses a parameter known as compactness which is equal to 
the number of PMTs in the bottom layer with more than 2 photoelectrons (PEs) divided by 
the number of PEs in the PMT with the largest number of PEs (in the bottom layer). A 
requirement that the compactness be greater than 2.5 retains 51% of the gamma rays while 
removing 91.5% of the cosmic ray background (for events with more than 50 PMTs in the 
trigger), resulting in an improvement in the gamma-ray flux sensitivity of 1.7. Compactness 
greater than 2.5 was required for this analysis. 

Milagro began data taking in 1999. The data set presented here begins on 15 December 
2000 and ends on 25 November 2003. Prior to 25 January 2002 the trigger was a simple 
multiplicity trigger requiring 60 PMTs to record a pulse greater than 1/6 of a PE within 180 
ns. After this date a risctime criterion was imposed on the cumulative timing distribution of 
struck PMTs. If a PMT signal is greater than 1/6 of a PE a trigger signal with an amplitude 
of 6.25 mV and width of 180 ns is generated. The signals from the 450 PMTs in the top layer 
are summed and sent to a VME trigger card. An 80 MHz flash analog-to-digital converter 
(FADC) digitizes the trigger signal. For each event that exceeds a "pre-trigger" condition 
the data from the FADC are stored in a FIFO. The risetime of the trigger signal is defined as 
the time taken for the trigger signal to go from 12.5% of its peak value to 88.5% of its peak 
value. For events with more than 75 PMTs no risetime requirement was imposed. Events 
with more than 52 PMTs were required to have a risetime less than 87.5 ns and events with 
greater than 26 PMTs were required to have a risetime less than 50 ns. The requirement 
on the risetime removes event triggers due to single muons at low multiplicity. The lower 
multiplicity requirement increases the effective area to low-energy events, which improves 
the sensitivity to gamma-ray bursts. The risetime requirement has a minimal impact on the 
sensitivity to Crab-hke sources. 

The trigger rate during this time varied from 1500 Hz to 1800 Hz. A total of 1009 days 
of data were acquired during this period for a detector "on-time" of 92%. The data are 
calibrated and reconstructed (to give the core position, direction, and information used in 
rejecting the background of the incoming particle) in real time. Except for selected regions 
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of the sky, only the reconstructed information is saved to disk. This analysis utilizes the 
reconstructed data set. The start date of the analysis is determined by the date that the com- 
pactness parameter and an improved shower-core fitter were included in the reconstructed 
data. 



3. Analysis Strategy 

Two maps of the sky are constructed: a signal map (comprised of the actual numbers 
of events coming from each bin in the sky) and a background map (comprised of an estimate 
of the cosmic-ray background from each bin in the sky). The maps are binned in 0.1 x 0.1 
degree bins. To estimate the background a technique called "direct integration" (Atkins ct al. 
2003) is used. The method makes use of the fact that the earth rotates and that the detector 
response is solely a function of local coordinates and time, and that the cosmic rays constitute 
an isotropic background. The underlying assumption of the method is that the shape of the 
detector response (in local coordinates) does not vary over the period during which the 
background is accumulated, which is two hours in this analysis. This method naturally 
accounts for rate variations in the detector and makes a high statistics measurement of the 
background (roughly 12 times as much background as signal is accumulated for each point 
in the sky). 

A bin of size 2.1 degrees in declination (5) by 2.1/ cos((5) degrees in right ascension is 
used to search the skymaps for evidence of a source of TeV gamma rays . This choice of 
bin size is based on the angular resolution of the detector as measured by the shadowing by 
the Moon of cosmic rays and Monte Carlo simulations of gamma-ray showers. The 0.1 x 0.1 
bins in the maps are summed to form these larger bins. Large bins are formed with centers 
at the center of each 0.1 x 0.1 degree bin by summing the contents of all the small bins that 
intersect the large bin. The summed data from the signal map and background map are then 
compared. The declination range of the search is between 1.1 degrees and 80 degrees. The 
lower limit is determined by the center of the first large bin with a low edge of zero degrees 
declination. The upper limit is determined by the rotation of the sky, above 80 degrees the 
hour angle interval used to generate the background is only 2.5 times the size of the signal 
region. To calculate the significance of each excess or deficit the prescription of Li and Ma 
(eq. 17) (Li and Ma 1983) is used. 
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4. Results of Survey 

Figure la shows the distribution of excesses (in standard deviations) for all bins in the 
map comprised of the entire 3 years of data. There is a statistically significant surplus of 
points in the sky with greater than 4 standard deviation excess. The bulk of this surplus can 
be attributed to two known sources of TeV gamma rays: the Crab Nebula and the active 
galaxy Mrk 421. Figure lb shows the distribution of excesses with 2-degree regions around 
the Crab Nebula and Mrk 421 removed. This distribution is consistent with expectations 
from random background fluctuations. A Gaussian fit to this distribution has a mean of 
-2.9x10"^ and a sigma of 0.987, consistent with the number of independent entries in the 
histogram. 

Figure 2 shows the map of the northern hemisphere in TeV gamma rays for the data set. 
The Crab Nebula (R.A.=83.64, Dec=22.01) and the active galaxy Mrk 421 (R.A.=166.11, 
Dec=38.21, z=0.03) are clearly visible in the map. The significance of the excess at the 
location of the Crab is 6.3 standard deviations and at the location Mrk 421 the significance 
of the excess is 4.4 standard deviations. (After correcting for instrumental dead-time and 
other known effects the results on the Crab Nebula result in a measured gamma ray rate of 
10.0 ± 1.4 events/day. This is slightly different from the result given in Atkins et al.(2003) 
of 10.7 ± 1.6 events/day due to the fact that different calibrations where used in the onhne 
reconstruction for some time periods.) Table 1 gives the location of all regions with an excess 
of greater than 4 standard deviations. Only the pixel with the largest significance from each 
independent region is listed in the table. Near the Crab Nebula the pixel with the greatest 
significance is at the location of the Crab Nebula, while near Mrk 421 it is 0.5 degrees 
from the position of Mrk 421. The number of these regions is consistent with the expected 
fiuctuations in the background given the large number of trials incurred in examining the 
entire sky. Therefore no claim is made that these regions are sources of TeV gamma rays, 
though they may be interesting regions for foUowup observations with the more sensitive 
ACTs. The Whipple Observatory has performed a follow-up observation of the region near 
R.A.=79.9 Dec=26.8 between November 2002 and January 2003 and has reported an upper 
hmit of 90 mCrab (Falcone et al. 2003), below the sensitivity of this survey. 

Both the energy response and the sensitivity to gamma ray sources of Milagro are 
dependent upon the declination of the source. Figure 3 shows the median energy of gamma 
rays that trigger Milagro (determined from Monte Carlo simulation) , satisfy the compactness 
criterion, and are reconstructed within 1.2 degrees of their true direction, averaged over a 
complete transit (i.e. 24 hours of observation) as a function of the declination of the source 
for several source spectral indices. The requirement that the direction of the particle be 
reconstructed within 1.2 degrees of its true direction is imposed to account for the bin size 
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used in the analysis (a 2.1 degree wide square bin has the same area as a 1.2 degree circular 
bin). 

Establishing upper limits to the gamma ray flux from any given point in the sky is 
straightforward for galactic (nearby) sources of TeV gamma rays. The prescription of Helene 
(Helene 1983) is used to calculate the confidence limits on the number of signal events from 
each region of the sky. Since the response of Milagro is energy dependent, the flux upper 
limits obtained from these data are dependent upon the energy spectra of the possible sources 
of TeV gamma rays. Figure 4 is a 2-dimensional map of the sky with the 95% C.L. upper 
limits to the flux given at each point. These limits are based on the assumption of source 
spectra proportional to £'^2.59 similar to the Crab Nebula (Aharonian et al. 2002; Atkins 
et al. 2003) at these energies). In order to translate the observed upper limits on the number 
of excess events from a given location into an upper limit on the flux of gamma rays, the 
detector response is normalized using the results from the Crab Nebula. This procedure 
accounts for the dead time of the detector, calibration errors, and other systematic effects. 
At declinations near 36 degrees (the latitude of Milagro) the 95% C.L. upper limits are on 
average 275 mCrab or F{>1 TeV)<4.8x 10~^^ cm~^ s~^. As the declination of the source 
increases or decreases from this value the source spends less time near zenith (where the 
atmospheric overburden is least and the response of Milagro is best) and the flux upper 
limits increase. At declinations of 5 degrees or 65 degrees the average 95% C.L. upper limits 
are of the order 600 mCrab {F{>1 TeV) < 1.05x 10~^^ cm~^ s~^). Figure 5 gives the factor 
by which these upper limits must be multiplied for potential sources with different spectral 
indices. The use of this flgure is best explained with an example. From Figure 4 obtain an 
upper limit from a location, R.A.=180., Dec=40 (F(>1 TeV)< 291 mCrab). To flnd the flux 
upper limit for a source with a differential spectral index of -2.0, find the curve for such a 
source in Figure 5 (the solid curve). At the declination of the source the y-axis value (0.41) 
is the amount by which the fiux upper limit from Figure 4 must be multiplied to give the 
fiux upper limit for this source (119 mCrab). 

For an extragalactic source one must also account for the absorption of the gamma 
rays due to interactions with the EBL (Primack et al. 2000; Stecker and De Jager 2002; 
Hartman et al. 1999). In the absence of a reliable measurement of the EBL, a model of 
the intensity and energy spectrum is used. Figure 6 shows the effect of the absorption due 
to the EBL on the upper limits given in Figure 4. In this figure the "baseline" model of 
Stecker and De Jager (Stecker and De Jager 2002) is used to calculate the effect of the 
absorption of energetic gamma rays and the source spectra are assumed to be proportional 
to i?"^-®^. Using a procedure similar to that described above for Figure 5, one can find the 
flux upper limit for sources at different redshifts. For example, using the same location as 
above (R.A.=180 Dec=40) but a source at a redshift of 0.03 (with a differential spectral 
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index of -2.59), the curve for z=0.03 in Figure 6 has a value of 2.4, resulting in a flux upper 
limit for this location and redshift of 700 mCrab. This upper limit is the normalization of 
the power law spectrum of the unabsorbed source at the top of the atmosphere. Assuming a 
source with an intrinsic spectrum (before absorption by the EBL) represented by, 

dE - 

then, for this source /q would be 2.4 times that of a galactic source with the same spectral 
index. (Note that to calculate the absolute luminosity of the source one must also multiply 
by the square of the distance to the source. This has not been accounted for in Figure 
6.) Because the absorption of TeV gamma rays by the EBL distorts the source spectrum 
before it reaches the earth one can not use Figure 6 in series with Figure 5. In general 
the effect on the sensitivity due to different source spectra is smaller for distant sources, 
since the EBL tends to make the distant sources look more alike regardless of their spectra. 
Figure 7 gives an example of the relative sensitivity to sources at different redshifts with 
different spectra. Figure 7a is similar to Figure 6 but the sources here are assumed to have a 
differential photon spectrum proportional to E^"^ and in Figure 7b a spectrum proportional 
to E~^. The figure gives the ratio of the flux upper limit for a source at the given redshift, 
declination, and spectrum to a local source (z=0.0) with a differential spectral index of -2.59. 
A perhaps surprising feature of this figure is that for more distant sources a source with a 
harder intrinsic spectrum is required to have a larger luminosity than a source with a softer 
spectrum for Milagro to make a detection. This is due to the effect of the EBL, where the 
high energy photons that the source is required to emit (by the model) are absorbed in 
transit and do not affect the ability of Milagro to observe the source but do count as part of 
the intrinsic source luminosity. 



5. Conclusions 

A complete survey of the northern hemisphere (dechnation 1.1° to 80°) for point sources 
of TeV gamma rays has been performed. These limits apply to the average flux level during 
the roughly 3 year period from December 15, 2000 through November 25, 2003. The average 
95% CL. upper limits range from 275 mCrab to 600 mCrab depending upon the declination 
of the source and arc over an order of magnitude more restrictive than previous limits. A 
prescription has been given to calculate the corresponding upper limits for sources with 
different spectra and for extragalactic sources. For sources with differential spectral indices 
of -2.0 the upper limits are 57% lower. For a source at a redshift of 0.03 the flux limits are 
a factor of 2.4 larger. While these hmits are the best available to date, Milagro has recently 
been completed with the construction of an array of 175 water tanks surrounding the central 
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reservoir. A comparable dataset, with this now complete Milagro detector, would improve 
these limits by a factor ~2. 
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ported by the National Science Foundation (under grants PHY-0070927, -0070933, -0075326, 
-0096256, -0097315, -0206656, -0302000, and ATM-0002744) the US Department of Energy 
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Fig. 1. — (a) The distribution of significances of the excesses and deficits in the analysis 
of the D.C skymap of the northern hemisphere, (b) The same data with 2 degree regions 
around the Crab Nebula and Mrk 421 removed. The dotted curves are the best fit Gaussians 
to the data. 




Fig. 2. — The northern hemisphere as seen in TeV gamma rays. At each point the excess is 
summed over a 2.1 degree by 2.1/ cos((5) bin and the significance of the excess in standard 
deviations is shown by the color scale. 
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Fig. 3. — The median energy of gamma-ray events tiiat trigger Milagro, pass the compactness 
cut, and are reconstructed with 1.2 degrees of their true direction as a function of source 
dechnation (this is equivalent to the 2.1 degree square bin used in the search). The response 
of Milagro is averaged over a complete source transit (i.e. one day's observation) and the 
source differential spectral index, ct, for each curve is given in the legend. 
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HA. 

Fig. 4. — The 95% CL. upper limits on the integral flux of gamma rays above 1 TeV (as- 
suming an _£;-2.59 flifferential photon spectrum) from each point in the northern hemisphere. 
The color scale on the right is in units of the flux from the Crab Nebula. To enhance the 
contrast of the flgure only declinations below 75 degrees are shown. 
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Fig. 5. — The effect of the differential spectral index on the upper hmits shown in Figure 
4. The y-axis gives the ratio of the flux upper limit for a source with a differential spectral 
index as indicated by the curve to a source with a differential spectral index of -2.59. The 
use of this figure is described in the text. 
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Fig. 6. — The effect of redshift on the upper hmits shown in Figure 4. These results assume 
a source spectrum proportional to The y-axis is the ratio of the flux upper limit for 

a source at the indicated redshift to a source at a redshift of zero. 
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Fig. 7. — The effect of redsfiift on tlie upper limits sfiown in Figure 4. Figure 7a assumes 
a source spectrum proportional to E~'^ and 7b a source spectrum proportional to E^^. In 
both cases the y-axis is the ratio of the flux upper limit for the described source (spectral 
index, redshift, and declination) and a local source with a differential spectral index of -2.59. 
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Table 1. The locations of all regions with greater than a 4cr excess. Only independent 
regions are entered in the table. Errors on the location of any possible source are ~0.5 
degrees. The two brightest points on the map are due to the two detected sources: the 
Crab Nebula^ and Mrk 421^. Upper hmits are not given for these two sources. The units of 
right ascension and declination are decimal degrees. The last column gives the 95% C.L. 
upper limit to the flux in units of the Crab flux. 
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